
    ##########################################################################################
    ### R code for Figure 12 of Supplemental Material
    ### Dynamic IRT Estimation Using CQ Key Vote Foreign Policy Votes
    ##########################################################################################

    library(MCMCpack)
    year.odd <- seq(from=1947, to=2010, by=2)
    T <- length(year.odd)

    Data <- read.table("Supplemental_Figure12_Data1.txt", header=T, sep="") # Combined CQ Key Votes (Foreign Policy ) Data - House
    Time <- read.table("Supplemental_Figure12_Data2.txt", header=T, sep="") # Time Index
    info2 <- read.csv("Supplemental_Figure10_Data3.csv", header=T) # Legislator Information: ICPSR_ID, congress, party ID, and Name (Adapted from Voteview.com data)

    out <- MCMCdynamicIRT1d(Data, 
                            item.time.map=Time$time,
                            mcmc=2000, burnin=500, thin=2,
                            verbose=1000, store.item=FALSE,        
                            theta.constraints=list(A14066A="+",    # Don Young (14066) 
                                                   A13035A="-"))   # Charles Rangel (D-NY)

    out_ideal <- out[,1:13407]  
    iter <- nrow(out_ideal)

    ### Prepare the MCMC chains to merge with legislator information
    original.label <- colnames(out_ideal)
    MC_Data <- data.frame(t(out_ideal))
    MC_Data$original.label <- original.label
   
    ### Decode the lable into ICPSR_ID and Congress
    congress <- icpsr_id <- time <- a22 <- a23 <- rep(NA, length(original.label))
    for(j in 1:length(original.label)){
    a22[j] <- strsplit(original.label, split="A.t")[[j]][1]
    time[j] <- as.numeric(strsplit(original.label, split="A.t")[[j]][2])
    icpsr_id[j] <- as.numeric(strsplit(a22, split="theta.A")[[j]][2])
    }
    congress <- time + 79
    
    Info_dynamic <- data.frame(original.label, icpsr_id, time, congress)
    colnames(Info_dynamic) <- c("original.label","ICPSR_ID","time","congress") 
    
    ### Merge with Legislator Information
    Info <- merge(Info_dynamic, info2, by=c("ICPSR_ID", "congress"), all=F)
    
    ### Merge Legislator Info and the Summary of Estimates
    Outputs <- merge(Info, MC_Data, by="original.label", all=F)
    info_columns <- 8 # the number of columns that are NOT MCMC chains

    ### Storages
    Dem.median <- matrix(NA, T, iter)
    Rep.median <- matrix(NA, T, iter)

    ### Loop over congresses

    for (t in 1:T){
    Est <- subset(Outputs, time==t)
    
    for (j in 1:iter){
    Dem.median[t,j] <- median(Est[,info_columns + j][Est$party==100])
    Rep.median[t,j] <- median(Est[,info_columns + j][Est$party==200])
    }
    }
    
    ### End of the loop
     
    Dem.median.quant <- apply(Dem.median, 1, FUN=quantile, na.rm=T, probs=c(0.025, 0.5, 0.975))
    Dem.median.mean <- apply(Dem.median, 1, FUN=mean, na.rm=T)
    Dem.median.CI <- rbind(Dem.median.mean, Dem.median.quant)
    
    Rep.median.quant <- apply(Rep.median, 1, FUN=quantile, na.rm=T, probs=c(0.025, 0.5, 0.975))
    Rep.median.mean <- apply(Rep.median, 1, FUN=mean, na.rm=T)
    Rep.median.CI <- rbind(Rep.median.mean, Rep.median.quant)
 
    ### Plot
    postscript("Supplemental_Figure12.eps", height=5, width=8, horizontal=F)
    par(mfrow=c(1, 1))
    plot(1:T, Dem.median.CI[1,], type='n', col='blue', ylim=c(-3,3), xaxt="n", xlab="", ylab="")
    points(1:T, Dem.median.CI[1,],  type='o', lty=1, col='blue', cex=.7)
    segments(1:T, Dem.median.CI[2,], 1:T, Dem.median.CI[4,], col='blue')
    points(1:T, Rep.median.CI[1,],  type='o', lty=2, col='red', cex=.7)
    segments(1:T, Rep.median.CI[2,], 1:T, Rep.median.CI[4,], col='red')
    axis(side=1, las=3, at=1:T, labels=year.odd, cex.axis=.8)
    legend('topleft', legend=c("Democrats", "Republicans"), lty=c(1, 2), col=c('blue', 'red'),cex=.8)
    dev.off()
       



